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"^ '. 

^\ ^ Motivated by a recent proposal of a Bethe approximation for the triangular Ising antiferromagnet 

0^ ' [Phys. Rev. B 56, 8241 (1997)], which seems to predict a disordered phase at any temperature in 

zero field, we analyze in some detail several mean-field like approximations for this model, namely 

the Bethe approximation itself, the cluster variation method and the hard-spin mean-field theory. 

^ ■ We show: (i) that the disordered phase predicted by the Bethe approximation is unphysical at 

low enough temperature because of a negative entropy; (ii) how the results of the cluster variation 
t^^ , method (namely, zero temperature entropy and critical temperature of the spurious transition) 

converge to the exact ones for increasing cluster size; (iii) that it is possible to construct a cluster 

(-H I variation approximation which yields a disordered phase which is stable down to zero temperature; 

Q_) . (iv) a few, so far unknown, zero temperature results (entropy and internal energy) of the hard-spin 

O ' mean-field theory. 
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PACS numbers: 05.50.-l-q, 64.60.Cn, 75.10.Hk 



I. INTRODUCTION 



The antiferromagnetic Ising model on the triangular lattice can be considered as the prototype of frustrated lattice 
'■^ models. It is defined by the (reduced) hamiltonian 

g: H/|JH^fT,a,-/i^a„ (1.1) 

■ (y) i 

^\\ ', where J < is the antiferromagnetic interaction, h = H/\J\ is the reduced uniform magnetic field and cTi is an Ising 
I> ' spin at site i, which can take the values -1-1 and —1. The former sum runs over nearest neighbor (NN) pairs and the 
C^ latter runs over all lattice sites. In zero field the ground state of the model is highly frustrated since the minimization 
^l I of bond energies would require that all pairs of NN spins are antiparallel, but this condition cannot be fulfilled even 
i^ ' on a simple triangle. On the contrary, six out of the eight possible configurations of a triangular plaquette (that is, 
f— ^ , all except the two fully parallel ones) are selected as the lowest energjt-Qpes, and this yields a huge degeneracy of the 
qn^ ■ global ground state. In zero field, the model has been solved exactlyld'B and the solution shows that the disordered 
^\ , phase is stable down to zero temperature and the zero temperature entropy per site (in the thermodynamic limit) is 
"^ S/N « 0.323066 x ks, where N is the number of sites and ks is the Boltzmann constant. 

2 Although an exact solution is available, a lot of work has focused on mean-field approximations, since the present 

Q . model can be viewed as a playground for testing approximate methods against exact results before applying them to 
more general and difficult models. The usual (single site) mean-field approximation fails qualitatively to predict the 
zero field behavior of the model: the disordered phase is thermodynamically stable only down to a certain transition 
temperature, bekiw which one finds a three-sublattice ordered phase. Also the more sophisticated Cluster Variation 
Q Method (CVM)ErQ gives the same qualitative picture, although with a lower transition temperature, which decreases 
upon increasing the cluster size. An interesting approach, which seeniKpto display a qualitatively correct behavior, 
is the Hard-Spin Mean-Field (HSMF) theory by Berker and coworkersB~l3, although quantities like entropy and free 
energy are to be calculated according to recipes which are not very well defincdj. Furthermore, a new investigation 
?H ' based on the Bethe-Pcicrls (BP) approximation has been recently carried out by Tamashiro and Salinastj: the 
authors claim that this approximation predicts a paramagnetic phase which in zero magnetic field is stable down to 
zero temperature. 

The purpose of the present paper is twofold: we shall first show that the paramagnetic phase found by Tamashiro and 
Salinasc2lis unphysical at low temperatures, by explicitly calculating the entropy, which turns out to be negative; then 
we shall address the more general issue of analyzing carefully the other approximations, with a particular attention 
to the zero field, zero temperature case. In the case of the CVM we shall show that the zero field, zero temperature 
entropy of the disordered phase converges, for increasing cluster size, to the exact value according to a power law, 
while the critical temperature converges very slowly to the exact value of zero; furthermore we shall shortly discuss a 
CVM approximation which does not predict a spurious transition in zero field. In the case of the HSMF theory we 



shall give a well-defined recipe to calculate the free energy and the entropy and use it to evaluate the zero field, zero 
temperature entropy. 

II. THE BETHE-PEIERLS APPROXIMATION 

The BP approximationliil has long been used to improve upon the ordinary mean-field results. It is the simplest 
mean-field like approximation which^;akes into account two-site correlations. 

Recently, Tamashiro and Salinasll3 applied the BP approximation to the present model, splitting the triangular 
lattice in the natural way into its three sublattices. They report that in the case of vanishing magnetic field the 
disordered (paramagnetic) phase is stable down to zero temperature. This is actually what the approximation predicts, 
but a deeper analysis is in order to understand what kind of disordered phase we are talking of at low temperature. 
Since in the next section we are going to discuss the CVM, it is worth mentioning that the CVM itself reduces, if one 
chooses a NN pair as the maximal cluster, to the BP approximation. Let us consider the reduced entropy density 
s — S/Nks given by the CVM for the triangular lattice-with the assumption of three non equivalent sublattices (say 
A, B, C) and the above choice of basic clusters (see e.g.ll3): 

s = - [Tt{pab InpAs) + Tr(pBC Inpsc) + Tr(pc'A InpcA)] 

+ ^ [Tt{pa \npA) + Tr(pB Inps) + Tr(pc Inpc)] • (2.1) 

Here Paa' is the NN pair density matrix (which is diagonal since we are dealing with a classical model) and its diagonal 
elements paa'icr^a') represent the probability of finding the NN spins belonging to the sublattices a and a' in the 
states a and a', respectively. Similarly pa = T^i-'a'iPaa') is the reduced (single site) density matrix, whose diagonal 
elements Pa{<7) represent the probability of finding a spin belonging to the sublattice a in the state a. Rewriting 
these probabilities in terms of the magnetizations ma = (ca) and the NN pair correlations Caa' = {caCTa') (where CTq 
and (7^ denote NN spin variables on the sublattices a and a', respectively and as usual (•) denotes thermal average) 
one obtains 



'1 + ania + (r'nia' + <Ta'caa'\ 5 V^ V^ ^f^ 

{aa') <y=±lcr'=±l ^ ^ a cr=±l 



2^ Z^ ^ Z +^2^2^^[ — — ' (2-2) 



where the outer sum runs over sublattice pairs {aa' = AB, BC\ CA) in the former term, over sublattices (a = A, B, C) 
in the latter and we have defined C{x) = x\nx. The (reduced) internal energy density (which more rigorously we 
shall refer to as an enthalpy) e ~ {Ti) /N\J\ can be written exactly as 

e= y^ Caa' - -y^^hmg. (2.3) 

{aa') a 

The minimization of the (reduced) free energy density 

g = G/N\J\^e~ts (2.4) 

(where t = fcsT/l J| and T is the absolute temperature) with respect to the magnetizations and the NN correlations 
can then be performed. We introduce the conditioned average of a spin (on sublattice a) with respect to its NN (on 
sublattice a') 

ma\a'{(^) ^ (ctq) = V a °" ' , (2.5) 

(Ja'=cr J^^ Pa'[cr') 

which is related to the previous variables by the identity: 

rriQ. + ctCqq' = TOQ,|Q/(cr)(l + ctTOq') cr = ±1. (2.6) 

The stationarity conditions for g can thus be written as 

^ = -l/.+ iuanh-V„ + iiVln^^^^m4 = (2-7) 

a'^a "1"^ ' 

dg tanh"^m„|a'(+)-tanh"^TOQ|Q/(-) 

-l + t - U, (_z!.8j 



dc, 



where the sum involves two sublattices only (for instance, ii a — A then a' — B,C). If the following quantity is 
defined 



A tanh ^ to„|q/(+) + tanh ^TOa|Q'(~) 

P'Hcta' ^ 



(being as usual /? — l/ksT), the equations |2.8| become 

■ma\a'{<^) = tanh (/3ryac</ -cr/3|J|) 



a = ±l. 



Substituting into Eq. 2.7 one obtains, with a bit of algebra: 

^ - rPH + -tanh^^mQ + ^ tanh"^ [tanh(/?| J|) tanh (/?77q,c,)] = 0. 



P\J\ 



dm a 



(2.9) 



(2.10) 



(2.11) 



a'T^Q 



This last equation, together with Eqs. 2.10, is equivalent to the self consistent equation derived in Ref. nO, and riaa' 
turn out to be the same effective fields as those of the BP approximation. Once the equivalence is proved we can go 
on working in the CVM approach. In the uniform case (ttIq = m Va and Ca^i = c Va 7^ a') the stationarity conditions 
for g can be easily written (without the above manipulations) in the following form: 



m = ±-\l{l + c)2 - (1 - c)2e-4/t 

_i 2m _i 

h = t 3 tanh 5 tanh m 

1 + c 



(2.12) 
(2.13) 



Eqs. 

be obtained 



Eq. 2.12 allows us to evaluate the magnetization m as a function of the pair correlation c and the reduced temperature t. 
This equation, in the limit t —^ 0, states that a disordered phase (tti = 0) can only have c = — 1. The zero temperature 
entropy and enthalpy (the latter coinciding with the internal energy u — U/N\J\ ii h = 0) can then be evaluated by 
I and 2.3 and turn out to be respectively s = —2 In 2 and (being h = 0)e — u = 3c= —3. The same results can 
by the BP free energy derived in Ref. |o[ If h — then the magnetization (and hence the local effective 
field) of the paramagnetic phase vanish, which corresponds to taking xi = X2 = in Eq. 24 of Ref. nfl. In this case 
the free energy reduces to (using v = tanh(/3J)) 



(i\J\g 



ln2 



■ln(l 



P- 



21n2-3/3|J|, 



(2.14) 



where the /3 ^ 00 asymptotic expression clearly agrees with our previous discussion. A couple of comments on these 
results are in order. A pair correlation equal to —1 means that every pair of NN spins is, in the ground state, in an 
antiparallel state. This is of course not possible, as mentioned in the Introduction. The BP approximation, taking 
into account only NN pairs and single sites, completely neglects frustration effects and predicts an intcrnaLenergy 
(and hence, at low enough temperatures, a free energy) which is much lower than the exact one, Uox — — lEJfl thus 
"stabilizing" the disordered phase. We note in passing that the BP approximation usually gives a free energy which 
is an upper bound to the exact one. 

It is clear that the state which is predicted to be the ground state cannot exist, and this is reflected in the negative 
value of the zero temperature entropy. Of course, such unphysical effects will not be limited to the point (t = 0, /i = 0), 
but will be found in a finite region of the phase diagram around this point. In order to give an idea of the size of 
this region we only mention that the zero field entropy vanishes at a temperature t w 1.28, and is negative below 
this temperature. Another, more stringent, criterion, can be based on the NN correlation c: a NN pair probability 
distribution can be obtained by partial trace from a triangle probability distribution only if c > — 1/3. In the BP 
approximation c = —1 in the ground state and c — —1/3 at a temperature t = 2/ln2 « 2.89. Below such a 
temperature the BP results are certainly unphysical, since they predict a NN pair probability distribution which 
cannot exist in a triangular lattice. In Fig. n we report the curves defined by s = and c = —1/3 in the (/i, t) p lane 
The former is evaluated numerically by Eq. 2.2 while the latter can be analytically drawn by substituting Eq. 
into Eq. |2.13|. 
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III. THE CLUSTER VARIATION METHOD 



In this section we develop a CVM approximation with the usual assumption of three non equivalent sublattices. 
We increase the cluster size anisotropically (just along one direction), according to Kikuchi's B2L-hierarchy, which has 



been shown to converge to the exact scJutionEJ. The proof of convergence has been made rigorous in the special case 
of the Ising model on a square latticecj. This approach is more convenient, from the point of view of computational 
efhciency, than choosing larger basic clusters in both directions, as in Refs. ||,Q. In our approximation the basic 
clusters are composed by L (downward pointing) triangles and contain 2L+ 1 lattice points, labelled as in Fig. y. Due 
to their shape we will denote such clusters as W-class. As usual the hamiltonian can be written as a sum of cluster 
hamiltonians: 

H/|J|=5]/.w(^l!Il.,2L+i), (3.1) 

w 

where the sum runs over all basic clusters and <t^ 2 2L+1 ^^ ^^^ spin configuration of the current cluster. ^From now 
on we shall use the following notation: 

""11,12, ■■■■in ^ i'-'^n ' '-^i2 1 ■ ■ • J '^in I ' V"^-^/ 

indicating a set of spins. The cluster hamiltonian has the following expression: 

_. 2L 1 2L-1 , 2L + 1 

h^ {^i,2,...,2L+i) = ^ E '^^^^+1 + 2(2L - 1) E "^"'^+^-2(21 + 1) E '^^ (3.3) 

i=l ^ ' i=\ ^ ' j=l 

where ui is the i-th spin (see Fig. 0), while the first sum runs over oblique bonds, the second over horizontal bonds 
and the third over sites. The coefficients before each sum avoid multiple countings: every oblique bond is shared by 
2L basic clusters, every horizontal bond by 2(2L — 1) and every site by 2(2L + 1) basic clusters. In the hypothesis 
of three non equivalent sublattices (labelled A, B and C) three non equivalent basic clusters (W^^\ W*^^) and W'^', 
where superscripts denote different spin probability distributions) must be taken into account, as shown in Fig. 
(first row). Observing that the same fraction (1/3) of each cluster type is present, the enthalpy can be written as 

^^sE E Pw^('^l,2,.^^,2L+l)^w(<Tl.2,^^^,2L + l), (3.4) 

^ — 1 ^l,2,...,2Z, + l 

where the k, index scans cluster types and p^ (<y\,2....,2L-\-i) is the probability of the spin configuration <Ti 2 2L+1 of 

the cluster W'^"-' (the inner sum runs over all possible configurations). Notice thatpthe coefficient 2 is the total number 
of W-clusters per lattice site. Following the most recent formulation of the CVM113 the entropy can be easily written 
as a linear combination (with suitable coefficients) of cluster entropies relative to a set of basic clusters (W^-'^^ , W^^^ , 
W^^-* in our case) and their subclustei-s. For a given cluster 7 the coefficient is obtained by the product of a number 
a.y, evaluated by Moebius inversiontj, and the ratio N-y/N (number of 7-clusters per lattice site). It turns out that 
in our problem only three subcluster classes have non vanishing a^. They are displayed in Fig. H (labelled by N, 0, 
E) while the corresponding coefhcients are given in Tab. |[ We finally obtain 



1 ' 



Y^ "^ „('') 



3 2^ 2^ PW ('^l,2,.^^,2L+l) [21n/}^^ ((Ti,2,...,2L + l) 



s = — 

-1 CTi 2,...,2L + 1 



Si^) 



-21np^''' {u\a,...ah)-^^^PQ (cri,3,^^^,2L+i) + lnpE ' (0-2^4, ■■■,2l) , (3.5) 



where the subcluster probability distributions have been defined: 



pV {crx^2,...aL) 


-Y.^^ 




C"2Z, + 1 


PO ((^1,3.....2L + 1) 


- >: p! 




^2.4 2i 


Pe \(^2A,....2L) 


- >; 




"■1,3 2i + l 



,2L+l) 



cri,2.,...,2L+i) 



Pw (0-1,2,. ■■,2L+l)- (3.6) 



It is then evident that the free energy g — e ^ ts is a function of the probability distributions of the basic clusters 
only. The variational procedure, with respect to these variables, has yet to satisfy three kinds of constraints: 



1. the normalization conditions (3 constraints) 

X! Pw' ('^i,2,...,2L+i) = 1 K=l,2,3; 



' 1,2,...,2L + 1 



2. the "translational" compatibiHty conditions (3 x 2^^ constraints) 

Pn ('^1,2,....2l) = Pn'' ('^1.2,...,2l) V(Ti_2,...,2L ; K= 1,2,3 

where 

PiM' (,'^2,3,...,2L + lJ — 2^Pw (,'^l,2,...,2L + lj 



(3.7) 



(3.8) 



(3.9) 



and the k indices are understood modulo 3 (as in following occurrences) 
3. the "rotational" compatibility conditions (3 x 2'' constraints) 



„(«) 



- >-i) 



where 



Pv ('^1:2,3) = Pv ('^3,1,2) Vcri_2,3 ; K = l,2,3 



Pv^ (0-1,2,3) = Yl Pw^ (cri,2,...,2L+l) . 

C4,5,...,2J: + 1 



(3.10) 



(3.11) 



The last two constraints can be more easily understood by the examples shown in Fig. y. Notice that a new subcluster 
class has been defined, namely V (the class of simple triangular plaquettes), while, due to constraints, N' is actually 
the same class as N . In order to take into account the above described constraints we have to introduce the Lagrange 
multipliers and define the free energy functional 



g = g- 



^Em^'^' E pIv^^ (-1,2....,2L + i) 



K=l 

3 



' 1,2,...,2Z, + 1 



E E -^^^k.2,...,2l) pi:H<ri.2,...,2L) ~ pi:r'H 



K=l Cri,2,-..,2L 



,""1,2,. ..,2L) 



k;=1 <7'i^2.3 



^ ^ u;(-^ (^1.2,3) p[;^ (Ti,2,3) - pr'^ (<T3,i,2) 



(3.12) 



where g is defined by Eqs. 2.4, 3.4 and [3.5| , while /i*-"-*, u''^-' {o'i.2,...,2l) and w^'^^ (''■i_2.3) are the Lagrange multipli- 
ers related to constraints 3.7, 3.8, 3.10 respectively. Taking the derivatives of g with respect to the probabilities 



Pw ('''i,2,...,2L+i) and setting them to zero, one easily obtains the Natural Iteration (NI) equationsi 



§ 



Pw (0-1,2,. ..,2L+i) = exp 
X exp 
X exp 



Ai(''V2-/iW (0-1,2,. ..,2L + l) 

W*"^ (0-1,2,. ..,2l) - W^''^^^ (0-2,3,. ..,2L+l) 

a;'"^(Tl,2,3)-Cc>(^+l)(<T2,3,l) 



PiM ' (0-1,2,. ..,2l) 



„(«) 



Pq iO-l,3,...,2L + lJ 
Pe (,0-2,4,...,2LJ 



-,1/2 



(3.13) 



where the "hat" on the left hand side denotes the evaluation of the probabilities at the "next step" of the NI method. 
The convergence is reached when a suitable distance between the current step and the next step turns out to be less 
than some tolerance. In order to compare results from different L, and hence with a different number of variables, we 
have used the 00-norm instead of the l-norm. The stopping test is thus 



max max 

K=l,2,3 (Ti2,-..,2L + 1 



~(«) 



Pw lcri,2,...,2L+lj - Pw (,cri,2,...,2L+lj 



M 



<e, 



(3.14) 



where usually e = 10 ^ has been chosen. The entire set of Lagrange multipliers can be determined at every step by 
means of a nested iterative procedureta. The general scheme of the method is as follows. 



A^^) 



1. Choose a set of guess values for p\^' ((Ti,2,...,2L+i) and set to zero the Lagrange multipliers v^'^' (<ti_2,...,2l) and 

uj^^Hf^ 1.2,3)- 



1-^-1 by means of Eqs. 3.12. 



2. Evaluate the unnormalized p{^ {ci.2,...,2L+i] 

/J, — u 

3. Compute a new estimate of the multipliers v and uj by means of the following equationsEj: 

Pn^ (cri,2....,2L) 



v^'''> {a- 



1,2,. ..,2L 



) = u(") {a- 



1,2 2L) 



b,, In 



PW \<^\.2,...,2L) 



cD^-) (cri,2,3) =u:^-^^{rT 



(«^) 



l,2,3j 



hu, In 



p(r^ (0-1,2.3) 
Pm ('^3,1.2) 



(3.15) 



where h^ and 6t^ must be empirically chosen in order to "stabilize" the numerical procedure (we have found that 
0.3 is a good value for both parameters). 

4. Repea t steps 2 and 3 until convergence is reached for the Lagrange multipliers (with a stopping test similar to 
(|3.14) but generally a different tolerance e'). 



5. Evaluate /i'-"' according to the normalization conditions (Eqs. |3.7[), which, applied to Eqs. 3.13, give 



X] Pw' (cri,2,...,2L+i; 

<Tl,2,...,2Z. + l 



fiC-J^O 



cxp 



(-M^'^Va) 



(3.16) 



and determine correctly normalized probabilities. 
6. Repeat steps 2-5 until convergence is reached for the probabilities, according to the stopping test ( |3.14 ) 
Finally, when convergence is reached, the equilibrium free energy can be easily evaluated as 

3 



9 



-sE"'-. 



(3J7) 



We have carried out the investigation increasing the cluster size parameter L, as far as convergence could be reached 
with a reasonable computational effort (L = 9 for the disordered phase and L = 8 for the ordered phases). It will be 
clear in the following that it is convenient to rewrite the sublattice magnetizations niA , rriB , nic as follows: 

mA = Re{m + Me*'^} 



ms = Re |to + Me'('^+2'^/^) | 



mc 



Re{TO + Me^('^+4'^/3)|, 



(3.18) 



where m and Me^^ are the new order parameters. The mapping can be inverted by 



raA + ras + mc 

771 = 



Me ^ = mA — m, + I- 



V3 



(3.19) 



where it is evident that m is the average magnetization and Me*'^ is an additional (complex) order parameter char- 
acterizing the breaking of translational symmetry. Several guess solutions have been tried: a uniform disordered one, 
with m and M equal to zero {ip undefined), and different symmetry-broken ones, with both zero and finite average 



magnetization m. The guess probability distributions have then been obtained in the hypothesis of statistically inde- 
pendent spins. First of all the zero field, zero temperature case has been considered. The disordered solution correctly 
displays magnetizations uia ~ rnB = rnc = 0, pair correlations cab — cbc — c-ca = ^1/3 and hence internal energy 
M = — 1 for any value of L (which is consistent with the fact that only the configurations having 2 spin "up" and 1 
spin "down" on each triangle or vice versa are allowed). Xhje entropy increases upon increasing L (as displayed in 
Tab. y) and converges to the exact zero temperature valueEla Sox ~ 0.323066. A power law behavior of the coherent 
anomaly type (see Ref. 17 and references therein) is observed: 

SL = s^-aL-^, (3.20) 

where sl is the entropy computed with a cluster size equal to L. Discarding the data points with L = 1 -;- 3 a least 
square fit to the above law has been obtained with Soo ~ 0.323126 (« Sex), « ~ 0.04647 and ■0 ~ 1.7512 (remarkably 
close to 7/4) , as shown in Fig. ||. 

As far as the ordered phases are concerned all guess solutions turn out to converge to three types only, which can 
be distinguished by the value of </?. To be more precise (p is exactly 7r/3 and 7r/6 (for any L) for two types of solutions 
and seems to approach 7r/4 (upon increasing L) for a third type. Taking into account the degeneracy induced by 
the symmetry of the hamiltonian we actually have ip ~ 7r/3 + fc7r/3 (with A; = 0, . . . , 5) or (/? = 7r/6 + fc7r/3 (with 
fc = 0, . . . , 5) for the former two types, and (f k, tt/A + kir/Q (with fc = 0, . . . , 11) for the third type. Particular care is 
needed in using the NI method for within the commonly used tolerances the three types of ordered solutions seem to 
be all present for any value of L but, if precision is increased, only one or two of the symmetry broken phases turn out 
to be real minima of the free energy. Numerical results are summarized in Tab. |l[ where the entropy and the order 
parameters of each phase (also the metastable ones) are displayed as a function of the cluster size L. The alternation 
between the different ordered solutions is evident. Also in this case the internal energy does not depend on L and turns 
out to be u = — 1, while a convergence of the entropy towards the exact value is recovered but no power law can be 
observed. On the contrary it is possible to see that the average magnetization rn of each ordered phase is always very 
small (actually vanishing for the phase with ip = 7r/6) and approaches zero upon increasing L (see Tab. |l|). The value 
of M turns out to decrease with L, too. The corresponding states in terms of sublattice magnetizations can be easily 
obtained by Eqs. 3. IS . It is also possible to observe (see Tab. |l|) that the difference of the entropies of the different 
minima gets lower and lower upon increasing the cluster size and this mimics the fact that the model has a continuous 
rotational symmetry (i.e. its free energy is independent of Lp), which is rclatcd_tp-.the existence of a field induced 
Kosterlitz-Thouless transition to a long range ordered phase in the ground statc£3'E£l. Notice that the entropy of the 
paramagnetic solution is always lower than that of the ordered solutions i.e. a (spurious) phase transition is predicted 
(at least for aay L investigated), but this is a common feature of CVM approximations (for instance the hexagon 
approximatioio predicts a stable phase of the type tp — 7r/6). Nevertheless we have verified that the approximate 
free energy tends to the exact one (which is known in zero field), upon increasing L. We have also calculated the 
critical temperature ic, limiting the analysis to L = 1 -^ 5, because computation is very time consuming near critical 
points. The transition temperature decreases upon increasing L, as expected, butj-the convergence towards the exact 
one (zero in this case) is not so fast as in other investigations based on the CVMlj. Anyway our results seem to be 
compatible with an asymptotic behavior of the following form (see Fig. H): 

6-2/*"= cx l/L (3.21) 

whepae"^'*, which in the low temperature limit is proportional to the correlation length, is the natural scaling variable 
hereEj. A more detailed investigation about the convergence of the CVM B2L approximation hierarchy towards the 
continuous symmetry present in the model and the convergence of the (spurious) transition temperature towards the 
exact one is beyond the scope of the present paper and is left for future work. 

Before concluding this section, devoted to the CVM, we shall briefly introduce a simple CVM-like approximation, 
known as "cactus" approximationEil, which turns out to be qualitatively correct for the Ising triangular antiferromag- 
net. In the cactus triangle approximation the entropy expansion takes into account only upward- (or downward-) 
pointing triangles as basic clusters and hence it can be written as 



/ , P123 (""1,2,3) 



Inpi23 (0-1^2,3) - -^^^^^^Pi (o'i) 



(3.22) 



where P123 (<Ti,2,3) is the triangle probability distribution, while pi {oi) is the site probability distribution, correspond- 
ing to the i-th sublattice. Notice that, due to the choice of one kind of triangles only (upward- or downward-pointing), 
there is no pair contribution to the CVM entropy. It is not difficult to prove that this kind of approximation 
only predicts a disordered phase down to zero temperature (in zero field) and its entropy turns out to be positive 
(s = ln(3/2) K, 0.4055). Both characteristics, even though obtained by a very simple (and actually not well justified) 
approximation, turn out to be qualitatively correct. 



IV. THE HARD-SPIN MEAN-FIELD THEORY 

In this section, after briefly reviewing the HSMF approximation, we perform on this basis a nearly analytical 
evaluation of the zero temperature entropy. This is a new result (exact within the HSMF theory) and, though lirmtcd 
to zero temperature, is an attempt to go beyond approximations like that proposed by Kabakgioglu and coworkeraJ for 
the evaluation of HSMF free energy. Incidentally it is also possible to obtain several informations about the model's 
behavior near zero temperature in this approximation. 

It is easy to show that the following exact relation holds: 

(^o) = /tanh ^'^^-i"- V (4.1) 

where ao is a spin variable at some site, and {(Tr}S^ ^rs its six nearest neighbors. As usual (•) denote thermal average. 
The above equation is one of the Callen identitiesc3 specialized for the antiferromagnetic model on a triangular lattice. 



The classical mean field approximation can be derived from Eq. 4.1 by replacing (tanh(-) ) by tanh(-). On the contrary 
the HSMF approximation consists in evaluating the, right hand side average in Eq. |4.l| by assuming the nearest 
neighbor spins {(Jr}r=i ^-re statistically independentQ. With this assumption it is easy to derive from Eq. 4.1 an 
equation for the magnetization of a homogeneous phase, which we are now interested in. The magnetization m is 
equal to the thermal average of any spin. We can then write 

m = (ar) ^P{crr^+l}~P{crr = -l} r = 0,l,...,6 (4.2) 

and hence 

1 i Tn 
PK = ±1} = ^— r-0,l,...,6 (4.3) 

where P{-} denotes the probability of the event described within curly braces. Let v'^ and i/~ be the number of 
nearest neighbor spins "up " and "down" respectively. In the hypothesis of statistical independence introduced above 



and making use of Eq. 4.3 it is possible to write the probability distribution of t^+ and i^ by the following formula: 



P{i^+ = n} = P{i^- = 6 - n} = p(f ) ( i4r^ ) n = 0,...,6 (4.4) 



where we have defined the binomial probability distribution 

,W(x)^(^)x"(l-x)^-". (4.5) 



Pn 

Noticing that ii v^ — n (and v^ — 6 — n) then 
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^ar = i^+ -ly' =2{n-3), (4.6) 



it is finally straightforward to rewrite Eq. 4.1 into 



- = Et-h^^^^^pi«)(i±^), (4.7) 



n=0 



which is the equation for m we were looking for. This equation can be numerically solved and the magnetization of 
the homogeneous phase can be computed at the given values of temperature t and field h. The equation is evidently 
invariant under the transformation {h ^> —h , m <-> —m) and we have then performed the calculation only for 
h > and for several temperatures: the results are displayed in Fig. |7|. Notice that we have only computed the 
magnetization of the homogeneous phase but we have not proved that it is the stable phase. So after calculating the 
entropy at i = and h — we will have to make sure that the model does not predict another (inhomogcneous) 
phase. The limit magnetization curve at t = 0, displayed in Fig. [7|, is a discontinuous (step) function. The plateau 
values of magnetization rrik {k = 1, 2, 3) can be rigorously computed observing that 



T ^ h — 2(n — 3) r, „ , „,, 

lim tann ^ = sgn [h — 2(n — 3)\, 



(4.8) 



where sgn is the sign function, which returns +1 , —1 or when its argument is positive, negative or zero respectively. 
By substituting the above equation into Eq. 4.7, for 2(fc — 1) < ft, < 2fc (fc = 1, 2, 3) one obtains 



fc+2 



(6) 



1 + m 



n—k-\-3 



(6) 



1 + m 



These are polynomial equations which can be written in the more compact form 

'1 + m' 



1 + "^ _ p(6) 
2 '=+2 



= 



by defining the cumulative distribution 



(4.9) 



(4.10) 



^r'(-) = E^-^(-) 



(4.11) 



and maki ng us e of Newton's binomial formula. It is now easy to find all the solutions of the polynomial (6-th degree) 
equations 4.10| by common numerical routines; mk is, out of the solutions of the k-th equation, the only one in the 
real interval [—1, 1]: 



mi « 0.1056 
m2 « 0.3213 
TO3 ~ 0.5562. 



(4.12) 



It is also easy to see that for ft > 6 Eq. 4.7 simply reduces to m = 1. The information about the magnetization 



plateau values is already sufficient for the zero temperature zero field internal energy w (t = 0, ?7i = 0) to be evaluated. 
The following thermodynamic identity must be used: 



u{t = 0, m = 0) = u(i = 0, TO = 1) - / ft(i = 0, TO)dTO. 

Jo 



(4.13) 



The integral can be easily evaluated for h{t — 0, m), implicitly defined via Eq. 4.7 in the linnt t — > 0, is a step function: 

ft(t = 0,TO) = 2fc Vm e (TOfc,TOfc+i), (4-14) 

where k = 0, . . . , 3 and TOq = —mi and m^ = 1. Moreover m (i = 0, to = 1) = 3 for if to, = 1 then the spin pair 
correlation is c = 1 and m = 3c in the triangular antiferromagnetic model. We finally obtain 



u{t = 0, TO = 0) = 2(mi + TO2 + TO3) - 3 « -1.0339. 



(4.15) 



It is a remarkable fact that the pair correlation is then c{t = 0, to = 0) ~ —0.3446, which is quite close to the exact 
value —1/3, but it does not respect the compatibility condition c > —1/3. 

In order to evaluate the zero temperature zero field entropy s (i = 0, to = 0) we need further manipulations. Again 
we will make use of a thermodynamic identity: 



s{t = 0, TO = 0) = s(i = 0, TO = 1) + / — (i = 0, m)dm, 

Jo c/f 



(4.16) 



where the entropy of the saturated system s (t — 0, m — 1) vanishes. As far as the evaluation of the integral is 
concerned we have to deriv e an explicit expression for the partial derivative of h{t,m) with respect to i in i = 0. To 
do that let us consider Eq. 4.7 and, for fc = 0, . . . , 3, solve it with respect to (ft — 2k) /t. With a bit of algebra we can 
write 



TO, - 



\^ tanh 



h{t,m) -2(n-3) (g^ 

. Pn 



h{t,m)-2k 



tanh 



n=0 
n^k+'i 



1 +TO 



(6) / 1 + m 

PkU 



(4.17) 



Let us now take the limit < ^ of each side of the above equation. Making use of Eq. 4.8 and 4.14 we easily obtain 
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lim 

t-.o 



h{t,m) — h{t = 0,to) 



fc+2 
n=0 



1 



= tanh 



Vm e (rnfc,mfe+i). 



n=k+A 



6) 



1 



(6) 
PkU 



(4.18) 



Notice that the left hand side of this equation is just the definition of the partial derivative we were looking for but on 
the right hand side a different expression is obtained for each interval (mfc, nik+i)- Expanding the inverse hyperbolic 
tangent and introducing definition (4.11), we finally have 



dh, „ , 1, 
^(^ = 0,n.)^-ln. 



1 



^rn _ (6) 

2 ^k+2 



1 



1 



\_rn _ (6) 
2 '=+3 



Vm e (TOfe,mfc+i). 



(4.19) 



The above result is graphically reported in Fig. ^ Due to singularities it is not easy to perform an integration of 
such a function by common nu meric al quadrature algorithms, but the form of the function permits the following 
manipulations. Substituting Eq. 4.19 into Eq. 4.16 and making the change of variable x = (1 + m)/2 it is possible to 
write 



s{t — 0,m 



3 r . 

o)-E / 

A;=0 '-"'^ 



fXk + l 



In 



Pt%i^) 



dx 



Xk + l 



In 



^i?3(-) 



dx 



(4.20) 



where Xk — [^ + m k)/2 for k ~ 1, . . . , 4 and xq = 1/2. It has been possible to introduce absolute values because only 
one solution of Eq. 4.10 lies in the physically meaningful region ttt, G [^1,1] and hence in that interval the following 
relation holds: 



1 



^rn _ (6) 

2 ^k+2 



1 



>0 



m > TOfc, 



where A; = 0, . . . , 4. The integrals of Eq. 4.2Q| can now be computed by means of the following formula: 
/ In |p(a::)|rfa; = Re < xlnag + /^(a^ — Zi) \hi{x — z^) — 1] > + const.. 



(4.21) 



(4.22) 



where p{x) is any (g-th degree) polynomial, Zi is its z-th (in general complex) root and ag is the highest power 
coefficient. The numerical problem is then reduced to the computation of polynomial roots. We have obtained 



s(t = 0,TO = 0) « 0.3869, 



(4.23) 



which is positive and not so far from the exact value. 

As a final check we have examined the possibility of symmetry-broken phases, as usual limiting our investigation to 
the case of a tripartite lattice. The unknowns of this problem are the magnetizations of the three sublattices, which 
we have called as usual m^, tib and mc- Along the lines of the derivation of the equation for m in the uniform case, 
and observing that if a site is of type A its nearest neighbors are three of type B and three of type C and the same 
for all permutations of A^ B, C, it is straightforward to write the following three equations: 



rua 



E E *^^^ ■ 



2{n' + n" - 3) (3) 

T Pn' 



1 



mb \ (3) 



1 



(4.24) 



with {a b c) — {A B C),{B C A),{C A B). A numerical solution of the system has shown that two symmetry-broken 
solutions with two equivalent sublattices (analogous to those obtained by Bethe approximation) exist in the vicinity 
of t = 0. As displayed in Fig. H, and in analogy to what happens for Bethe approximation, symmetry-broken solutions 
are not present at ft- = because they become critical in two simmetrical points h/t — ±^. This result, as well as 



10 



the existence of two points {h/t = ±x) in which the uniform solution crosses one simmetry-broken solution, has been 
extensively p ointe d out in Ref. |^. Notice that in Fig. ^ magnetizations are functions of h/t only. This is due to the 
form of Eqs. |4.24| in which, for i -^ and in a region of h close to zero, only one hyperbolic tangent (the one whose 
argument is h/t) is significantly different from ±1. In this way we can be sure that the one displayed in Fig. y is the 
asymptotic behavior of magnetizations for i — > and hence no symmetry-broken solution can appear at /i = for any 
arbitrarily small temperature value. 

V. CONCLUSIONS 

We have investigated several mean-field like approximations for the antiferromagnetic Ising model on the triangular 
lattice. We have shown that the paramagnetic phase predicted by the BP approximation is unphysical in a low 
temperature, low field region since, due to the neglect of frustration effects, the approximation predicts as the ground 
state a state which cannot exist on the triangular lattice, and hence has negative entropy. In the case of the CVM 
we have investigated some of the convergence properties of the B2L series of approximations, showing that the zero 
field zero temperature entropy of the (metastable) disordered phase converges to the exact value with a power law of 
the coherent anomaly type, sl = Sox — aL^^ , with ip « 7/4. The CVM predicts a variety of ordered ground states 
(and hence a spurious phase transition) whose symmetry depends on the cluster size and a possible interpretation of 
this behavior has been suggested. Subsequently the critical temperature has been evaluated, limiting the analysis to 
L = 5, due to the well-known computational difficulties near critical points. The transition temperature turns out to 
decrease with L, as expected, although very slowly. We have also mentioned that the cactus triangle approximation 
behaves in a qualitatively correct way. Finally, we have considered the HSMF theory, calculating for the first time 
the zero field zero temperature entropy and internal energy by means of a well-defined procedure. 
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N 


-1 


2 





-1 


1 


E 


+1 


1 



TABLE II. Results obtained by the CVM for different cluster sizes (first column). The second and third columns contain 
respectively the argument <f (discriminating the type of symmetry breaking, ^p is undefined for the disordered phase) and the 
modulus M of the complex order parameter, whereas the fourth column contains the average magnetization m. The fifth 
column contains the (reduced) entropy s: for each value of the cluster size L the different solutions are ordered on decreasing 
entropy (the first is the stable one). For L = 9 only the disordered phase has been investigated. 



L 


ip (dcg.) 


M m 


s 


1 


30.0000 
60.0000 


0.89810 0.00000 
0.86914 0.00533 
0.00000 0.00000 


0.3095711 
0.3051012 
0.2876821 


2 


60.0000 


0.79214 0.00252 
0.00000 0.00000 


0.3177109 
0.3095711 


3 


30.0000 


0.70585 0.00000 
0.00000 0.00000 


0.3207863 
0.3164601 


4 


47.2021 
60.0000 


0.68361 0.00041 
0.67748 0.00045 
0.00000 0.00000 


0.3210938 
0.3210694 
0.3190285 


5 


60.0000 
30.0000 


0.65295 0.00046 
0.65240 0.00000 
0.00000 0.00000 


0.3217165 
0.3217141 
0.3203469 


6 


45.1102 


0.62727 0.00023 
0.00000 0.00000 


0.3221726 
0.3211130 


7 


30.0000 
60.0000 


0.59927 0.00000 
0.59818 0.00018 
0.00000 0.00000 


0.3223160 
0.3223155 
0.3215873 


8 
8 


60.0000 
30.0000 


0.58114 0.00016 
0.58017 0.00000 
0.00000 0.00000 


0.3224757 
0.3224749 
0.3219084 


9 


— 


0.00000 0.00000 


0.3221350 
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FIG. 1. In the field-temperature {h,t) plane the loci of points in which the BP approximation predicts s — (thick dashed 
line) and c = —1/3 (thick solid line) are displayed. Some contour lines of the NN pair correlation c are also reported (thin 
lines) . 

FIG. 2. Basic clusters of the B2L hierarchy on a triangular lattice. Numerical ordering of sites is that used in formulas. 

FIG. 3. Cluster involved in the B2L entropy expansion for a triangular lattice split into 3 non equivalent sublattices {A, 
B, C). Due to lattice splitting each cluster class (W, N, 0, E) is split into 3 subclasses or "types" { \ '■^•', '■^') with the same 
shapes but different probability distributions. 

FIG. 4. Examples of compatibility conditions, (a) For L = 2 it is shown that the cluster obtained by the superposition 
of W'^) {BCABC) and W*^' {ABCAB) is N'^' [BCAB, thick hues); hence the same probability distribution p^^' must be 
obtained either by tracing p^ over site 5 (C) or p^ over site 1 (A), (b) Again for L = 2 it is shown that the cluster obtained 
by the superposition of W'^^ (BCABC) and a rotated W'^) (ABCAB) is V'^' (BCA, thick lines), coinciding with V^^' (ABC); 
hence the same probability distribution p^ must be obtained either by tracing p^ over sites 4,5 {B,C) or p\fj over sites 
4,5 (A,B). It is easy to generalize both examples to any cluster size L. The same constraints must be satisfied by any pair 
(W'^^'W'^', W"^W''^', W'-^^W"^): conditions (a) imply a translational invariance and conditions (b) a rotational invariance. 

FIG. 5. The difference between the extrapolated zero temperature zero field entropy Soo ~ Scx and that calculated in 
the CVM approximation sl (for the disordered phase) is plotted versus the reciprocal cluster size 1/L. Dots represent data 
obtained by the CVM while the solid straight line is the least square fitting (data point L = 1 4- 3 discarded) . The slope of the 
straight line is the exponent tp. 

FIG. 6. The scaled critical temperature exp( — 2/ic) is plotted versus the reciprocal cluster size 1/L. Dots represent data 
obtained by the CVM while the solid straight line is a possible asymptotic regime converging to tc — 0, obtained as the tangent 
(in 1/L = 0) of a parabolic fit (thin dashed line). 

FIG. 7. Magnetization m vs. field h for the paramagnetic phase at different temperatures t = 0.1,0.4,0.7, 1.0 according to 
the HSMF th eory. The limit behavior t = is also displayed. Numerical values of the plateau magnetizations mi, m2, rris are 



given by Eqs. 4.1 J 



FIG. 8. Partial derivative of the field h with respect to temperature i in i = as a function of magnetization m. In analogy 
with Fi g. |7t m is reported on the vertical axis (and dh/dt on the horizontal axis). Numerical values of nii, 1x12, m-j, are given 



by Eqs. 4.15 



FIG. 9. Magnetization m, as a function of h/t in the asymptotic regime i — > (in the neighborhood of /i = 0). For both 
symmetry-broken solutions the magnetizations of the two non-equivalent sublattices are denoted by tub ~ rnc and tua', the 
two solutions become critical at h/t = ^. The uniform magnetization (solid line) crosses that of one broken-symmetry solution 
(dashed lines) at h/t — x- The other broken-symmetry solution is represented by dash-dotted lines. 
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